For execution on a local, multicore CPU with excess RAM we recommend calling
options(mc.cores = parallel::detectCores()).
To avoid recompilation of unchanged Stan programs, we recommend calling
rstan_options(auto_write = TRUE)
For within-chain threading using `reduce_sum()` or `map_rect()` Stan functions,
change `threads_per_chain` option:
rstan_options(threads_per_chain = 1)
Loading required package: gdata
Attaching package: 'gdata'
The following object is masked from 'package:stats':
nobs
The following object is masked from 'package:utils':
object.size
The following object is masked from 'package:base':
startsWith
Loading required package: bayesplot
This is bayesplot version 1.10.0
- Online documentation and vignettes at mc-stan.org/bayesplot
- bayesplot theme set to bayesplot::theme_default()
* Does _not_ affect other ggplot2 plots
* See ?bayesplot_theme_set for details on theme setting
Loading required package: stringr
Loading required package: dplyr
Attaching package: 'dplyr'
The following objects are masked from 'package:gdata':
combine, first, last
The following objects are masked from 'package:stats':
filter, lag
The following objects are masked from 'package:base':
intersect, setdiff, setequal, union
Loading required package: ggplot2
Loading required package: hBayesDM
Warning in library(package, lib.loc = lib.loc, character.only = TRUE,
logical.return = TRUE, : there is no package called 'hBayesDM'
datapath <-'/fast/work/groups/ag_schlagenhauf/B01_FP1_WP2/ILT_Stan_Modeling'out_path <-'/fast/work/groups/ag_schlagenhauf/B01_FP1_WP2/ILT_Stan_Modeling/Output'behavpath <-'/fast/work/groups/ag_schlagenhauf/B01_FP1_WP2/ILT_DATA'# load files containing true parameters used as input for simulationorig_file <-'fit_n_142_2024-01-29_bandit2arm_delta_PH_withC_estimation1_delta0.9_stepsize0.1.rds'orig_fit <-readRDS(file.path(out_path, orig_file)) # Stan model output# load file containing true parameters computed as transformed parameters during simulationsim_file <-'sim_2024-01-31_bandit2arm_delta_PH_withC_sim.rds'sim_fit <-readRDS(file.path(out_path, 'Parameter_Recovery', sim_file))# get file names for simulated data resultsrecovery_file <-'recovery_2024-02-01_bandit2arm_delta_PH_withC.rds'recovery_fit <-readRDS(file.path(out_path, 'Parameter_Recovery', recovery_file)) # Stan model outputcolor_scheme_set("mix-blue-pink")
2 Posterior predictive checks & correlations
# Load true parameters## extract posterior means for all parameters to use them as input for simulation### posterior means of parameters as input for simulationtrue_mu_pr <-as.vector(summary(orig_fit, pars="mu_pr")$summary[, c("mean")]) true_sigma <-as.vector(summary(orig_fit, pars="sigma")$summary[, c("mean")]) true_A_pr <-as.vector(summary(orig_fit, pars="A_pr")$summary[, c("mean")]) true_tau_pr <-as.vector(summary(orig_fit, pars="tau_pr")$summary[, c("mean")]) true_gamma_pr <-as.vector(summary(orig_fit, pars="gamma_pr")$summary[, c("mean")]) true_C_pr <-as.vector(summary(orig_fit, pars="C_pr")$summary[, c("mean")]) ### transformed parameters, output of simulationtrue_sim_pars <-as.vector(extract(sim_fit))true_A <-as.vector(true_sim_pars$A)true_tau <-as.vector(true_sim_pars$tau)true_gamma <-as.vector(true_sim_pars$gamma)true_C <-as.vector(true_sim_pars$C)true_mu_A <-as.vector(true_sim_pars$mu_A)true_mu_tau <-as.vector(true_sim_pars$mu_tau)true_mu_gamma <-as.vector(true_sim_pars$mu_gamma)true_mu_C <-as.vector(true_sim_pars$mu_C)## extract parameter values based on simulated datarecovered_mu_pr <-as.matrix(recovery_fit, pars ="mu_pr")recovered_sigma <-as.matrix(recovery_fit, pars ="sigma")recovered_A_pr <-as.matrix(recovery_fit, pars ="A_pr")recovered_tau_pr <-as.matrix(recovery_fit, pars ="tau_pr")recovered_gamma_pr <-as.matrix(recovery_fit, pars ="gamma_pr")recovered_C_pr <-as.matrix(recovery_fit, pars ="C_pr")recovered_A <-as.matrix(recovery_fit, pars ="A")recovered_tau <-as.matrix(recovery_fit, pars ="tau")recovered_gamma <-as.matrix(recovery_fit, pars ="gamma")recovered_C <-as.matrix(recovery_fit, pars ="C")recovered_mu_A <-as.matrix(recovery_fit, pars ="mu_A")recovered_mu_tau <-as.matrix(recovery_fit, pars ="mu_tau")recovered_mu_gamma <-as.matrix(recovery_fit, pars ="mu_gamma")recovered_mu_C <-as.matrix(recovery_fit, pars ="mu_C")